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Introduction tc 
Mission Desiai 


anetary 


■ The interplanetary design problem is composed of both discrete and 
real-valued decision parameters: 


Choice of destination(s), number of planetary flybys, identities of flyby 
planets 


Launch date, flight time(s), epochs of maneuvers, control history, flyby 
altitudes, etc. 


■ For example, for a main-belt asteroid mission, the designer must 
choose: 


- The optimal asteroid from a set of scientifically interesting bodies 
provided by the customer 

- Whether or not to perform planetary flybys on the way to the main belt 
and, if so, at which planets 

- Optimal trajectory from the Earth to the chosen asteroid by way of the 
chosen flyby planets 
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Unique Chara< 
Interplanetary 


■ Low-thrust electric propulsion is characterized by high power requirements but 
also very high specific impulse (l sp ), leading to very good mass fractions 

■ Low-thrust trajectory design is a very different process from chemical trajectory 
design 

- Chemical thrusters fire for minutes during a mission time-scale of years 

■ Chemical maneuvers may be approximated as occurring instantaneously (impulsively) 

■ This allows a mission designer to parameterize the cost of a chemical mission in units of “change 
in velocity,” or Av 

■ Av is invariant to changes in spacecraft hardware or launch vehicle, and mass fraction may be 
computed using only l sp and Av 

Electric thrusters fire for months or years, sometimes the entire mission duration 

■ The mission designer therefore must choose a time-history of thrust control rather than discrete 
maneuvers 

■ The impulsive approximation does not apply, and so the affect of the thruster must be continuously 
integrated as the spacecraft flies 

■ The propellant flow rate required to provide a given thrust is dependent on the available power at 
that instant 

"Therefore different choices of power and thruster systems change the trajectory solution! 
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Traditional M 
Flyby Traject 


rust, Multi- 


■ Several methods of picking the destination and flyby sequence: 

- Grid search over all possible choices of destinations, flyby sequence, 
propulsion system, power system, etc. (very expensive and often 
impractical) 

- Intuition-guided manual design of the trajectory (even more expensive, 
can miss non-intuitive solutions) 

■ Several methods of designing the trajectory: 

- Local optimization from an initial guess provided by a chemical 
mission design (but sometimes the optimal chemical trajectory does 
not resemble the optimal low-thrust trajectory) 

- Local optimization from an initial guess provided by a low-fidelity 
approximation to the low-thrust model, i.e. shaped-based methods 
(but sometimes the shape-based method cannot accurately 
approximate the true trajectory) 
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lanetary 


■ Gage, Braun, and Kroo, 1994 - autonomous chemical design with variable mission 
sequence (no deep-space maneuvers) 

■ Vasile and de Pascale, 2005 - autonomous chemical design for fixed mission sequence 

■ Vinko and Izzo, 2008 - autonomous chemical design for fixed mission sequence 

■ Wall and Conway, 2009 - autonomous low-thrust design for fixed mission sequence (no 
planetary flybys) 

■ Chilan and Conway, 2009 - autonomous low-thrust and chemical design for fixed mission 
sequence (no planetary flybys) 

■ Yam, di Lorenzo, and Izzo, 201 1 - autonomous low-thrust design for fixed mission 
sequence 

■ Abdelkhalik and Gad, 201 1 , 2012, and 2013 - autonomous chemical design with variable 
mission sequence 

■ Englander, Conway, and Williams, 2012 - autonomous chemical design with variable 
mission sequence 

■ Englander (dissertation) 2013 - autonomous low-thrust design with variable mission 


sequence 
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Automated Mi! 
Optimal Contr 


brid 


■ Break the mission design problem into two stages, or “loops” 

- “outer-loop” picks sets of destinations, planetary flybys, sizes the 
power system, can pick propulsion system - a discrete optimization 
problem 

- “inner-loop” finds the optimal trajectory for a given candidate outer- 
loop solution - a real-valued optimization problem 

- For the outer-loop to work, the inner-loop must function autonomously 
(i.e. no human interaction) 
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Multi-Objecti 


Control 


■ The customer (scientist or project manager) most often does not want 
just one point solution to the mission design problem 

■ Instead, an exploration of a multi-objective trade space is required 

■ For a typical main-belt asteroid mission the customer might wish to see 
the trade-space of: 

- Launch date vs 

- Flight time vs 

- Deliverable mass 

- While varying the destination asteroid, planetary flybys, solar array 
size, etc 

■To address this question we use a multi-objective discrete outer-loop 
which defines many single objective real-valued inner-loop problems 
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Outer-Loop T 


Optimization 


■ The outer-loop finds the non-dominated trade surface between any set 
of objective functions chosen by the user 

■ Non-dominated surface means “no point on the surface is superior to 
any other point on the surface in all of the objective functions” 

■ The outer-loop solver may choose from a menu of options for each 
decision variable 

■ The choices made by the outer-loop solver are used to define trajectory 
optimization problems to be solved by the inner-loop 
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Anatomy of 


• Break mission into a set of “journeys,” 
each of which in turn is broken into 
“phases” 

• The endpoints of a journey are chosen in 
the problem assumptions 

• The endpoints of a phase (i.e. a flyby 
target) may be chosen by the user or 
Outer-Loop solver 
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Outer-Loop ' 

An Examnle 

1 ■ 1 1 ha — . A III wj 1 

Asteroid Choices 


Flight time upper bound 



Array Size 


Flyby Choices (Journey 1) 


Flyby Choices (Journey 2) 


Code 

Body 


Code Flight Time (y) 


Code 

Power (kW) 


Code 

Body 


Code 

Body 


0 

Ceres 


0 

6 


0 

10 


0 

Earth 


0 

Mars 


1 

Pallas 


1 

6.5 


1 

11 


1 

Mars 


1 

none 


2 

Juno 


2 

7 


2 

12 


2 

none 





3 

Vesta 


3 

7.5 


3 

13 


3 

none 





4 

Astraea 


4 

8 


4 

14 








5 

Hebe 


5 

8.5 


5 

15 








6 

Iris 


6 

9 


6 

16 








7 

Flora 


7 

9.5 


7 

17 








8 

Metis 


8 

10 


8 

18 








9 

Hygiea 


9 

10.5 


9 

19 








10 

Parthenope 


10 

11 


10 

20 








11 

Victoria 


11 

11.5 











12 

Egeria 


12 

12 












13 

14 


Irene 

Eunomia 


15 

16 


Psyche 

Thetis 



Array Size 

Flight Time Upper 

Sample Mission 

Asteroid 1 Potential Planetary 

Asteroid 2 

Potential 

Code 

8 

Bound 

4 

0 

Flyby 1 

1 

1 

Planetary Flyby 2 

1 

Translation 

18 kW 

By 

Ceres 

Mars 

Pallas 

none 
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Multi-Objec 


■ The outer-loop optimization problem is solved using a discrete multi- 
objective solver, in this case Non-Dominated Sorting Genetic Algorithm 
II (NSGA-II) 

■ NSGA-II finds the non-dominated front, surface, or hyper-surface 
between any number of objectives chosen by the user 


Initial generation 









Population evolves via 
genetic operators 

> 


Mass 











TOF 


TOF 
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Inner-Loop M 


ation 


■ The inner-loop solves a real-valued trajectory optimization problem 
which is defined by each candidate solution to the outer-loop problem 

■ The inner-loop must function autonomously because the problems are 
generated in real time and there is no opportunity for human 
intervention 

■ The outer-loop is only as good as the solutions to the inner-loop 
problem, so the inner-loop must be robust 

■ A given run of the outer-loop may require hundreds or even thousands 
of runs of the inner-loop, so the inner-loop must be fast 

■ If the individual inner-loop runs are independent then many of them can 
be run in parallel 
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Multiple Gravit 
(MGALT) via th 


w-Thrust 
n Transcription 


■ Break mission into phases. Each phase starts and ends at a body. 

■ Sims-Flanagan Transcription 

- Break phases into time steps 

- Insert a small impulse in the center of each 
time step, with bounded magnitude 

- Optimizer Chooses: 

■ Launch date 

- For each phase: 

■ Initial velocity vector 

■ Flight time 

■ Thrust-impulse vector at each time step 

■ Mass at the end of the phase 

■ Terminal velocity vector 

■ Assume two-body force model; propagate by solving Kepler’s problem 

■ Propagate forward and backward from phase endpoints to a “match point” 

■ Enforce nonlinear state continuity constraints at match point 

■ Enforce nonlinear velocity magnitude and altitude constraints at flyby 



O patch point 
> burn point 
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Power, Propu 
Modeling 


■ Medium-fidelity mission design requires accurate hardware modeling 

■ Launch vehicles are modeled using a polynomial fit 

m delivered = (l — G Lv) ( a LV^I + ^LV^3 + C LV ^ 3 + d LV C 2 + e LV ^ + hv) 
where o LV is launch vehicle margin and C 3 is hyperbolic excess velocity 

■ Thrusters are modeled using either a polynomial fit to published thrust and mass flow rate data 

7 h= a F P 4 + b F P 3 + c F P 2 + d F P + e F 
T — QpP^ + bpP 2 + CpP 2 + d T P + 6 f 


or, when detailed performance data is unavailable 

2 TjP 


T = 


Isp9o 


Power is modeled by a standard polynomial model 


ro + y + %7 , 

1 (1 - T)‘ 


r 2 \ 1 + y 3 r + y 4 r 2 

where P 0 is the power at beginning of life at 1 AU and r is the solar array degradation constant 


Ephemeris data for solar system bodies is provided via the SPICE toolkit 
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Inner-Loop S 
Nonlinear Pr< 


Minimize f(x) 

Subject to: 

%lb — % — %ub 

c(x ) < 0 

Ax <0 

where: 

x w , x ub are lower and upper bounds on the decision variables 
c(x) is a vector of nonlinear constraints 
Ax is a vector of linear constraints 

■ There are several third party solvers that do this (SNOPT, IPOPT, fmincon, 
vfl 3AD) 

■ But all of these methods require an initial guess... 
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Inner-Loop S 
Monotonic B; 



Leary, 2000 

Vasile, Minisci, and Locatelli, 2009 
Yam, di Lorenzo, and Izzo, 2011 
Englander (dissertation), 2013 
Casioli et al., 2013 
Englander and Englander, 2014 


Improved from standard MBH by: 

1 . “Feasible point finder” aggregate penalty method 

2. Non-uniform (Pareto) perturbation step 

3. “Time-hop” operator (Casioli et al.) 
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Example: Ma 


Mission Objective 

Visit two main-belt asteroids with diameter greater than 50 km 
(475 bodies meet this filter) 

Launch Vehicle 

Atlas V 401 

Power System 

Array power at 1 AU 
Cell performance model 

15 kW 
1/r 2 

Spacecraft bus power 

800 W 

Power margin 

15% 

Propulsion System 

Thruster 

NEXT (throttle table 11, high-lsp mode) 

Number of thrusters 

1 

Duty cycle 
Propellant tank 

90% 

unconstrained 

Mission Sequence 

up to two planetary flybys are permitted before the first asteroid and up to 
one between the first and second asteroids 

Inner-Loop Objective Function 

Maximize delivered mass to second asteroid 

Outer-Loop Objective Functions 

Delivered mass to second asteroid 
Launch year 
Flight time 
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Launch Year 

Code 

Year 

0 

2020 

1 

2021 

2 

2022 

3 

2023 

4 

2024 

6 

2025 

7 

2026 

8 

2027 

9 

2028 

10 

2029 


Second Asteroid 

Code 

Body 

0 

Ceres 

1 

Pallas 

2 

Juno 

3 

Vesta 

4 

Astraea 

5 

Hebe 

6 

Iris 

7 

Flora 


(475 

choices) 


First Asteroid 

Code 

Body 

0 

Ceres 

1 

Pallas 

2 

Juno 

3 

Vesta 

4 

Astraea 

5 

Hebe 

6 

Iris 

7 

Flora 


(475 

choices) 


Flight Time Upper 


Bound 

Code 

# Years 

0 

5 

1 

6 

2 

7 

3 

8 

4 

9 

5 

10 

7 

11 

8 

12 


First Journey First Flyby 

Code 

Body 

0 

Earth 

1 

Mars 

2 

Jupiter 

3 

No flyby 

4 

No flyby 

5 

No flyby 


First Journey Second Flyby 

Code 

Body 

0 

Earth 

1 

Mars 

2 

Jupiter 

3 

No flyby 

4 

No flyby 

5 

No flyby 


Second Journey Flyby 

Code 

Body 

0 

Earth 

1 

Mars 

2 

Jupiter 

3 

No flyby 

4 

No flyby 

5 

No flyby 


1 . 1 6x1 0 9 
possible 
combinations, 
4.82x10 s with 
duplicates 
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Delivered mass 


Main-Belt Twc 
First Generati 
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Main-Belt Twc 
Final Generati 
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Main-Belt Tv 
T rajectories 



unpowered flyby 
Ea rth 
4/2/2021 

= 5,167 km/b 
DEC = 8.9* 
altitude = 300 km 

m = 2808 kg 


launch 

Earth 

2/6/2020 

C3 = 0.915 km 2 /? 
DLA - 13. 3 C 
m = 2976 kg 


LI rendezvous 
21Lutetia 
2/6/2032 
m = 1935 kg 


IT rendezvous 
270Anahita 
4/12/2028 
m = 2093 kg 


le8 


x (km) 

i- 


1 eI 
0 ^ 
-1 
-2 
-3 




departure 
270Anahita 
4/12/2029 
m = 2093 kg 


leS 


x (km) 


-^3 I 2 “ . I ? } l l t 



LT rendezvous 
172Baucis 
10/9/2029 
m = 837 kg 


LT rendezvous 
63Au Sonia 
11/26/2027 
m =* 901 kg 


departure 
63Ausonia 
11/25/2028 
m = 901 kg 


-2 

-3 

< 

-4 f 

\ 


launch 

Earth 

10/9/2024 
C3 = 33.762 km- /s 2 
DLA = -28,5° 
m = 1236 kg 


A 12-year mission to Anahita and 
Lutetia delivers a very large science 
payload 


A 5-year mission to Ausonia and 
Baucis delivers a smaller payload in 
less time and at a lower cost 


Both of these candidate missions, and many others, would be 
valuable data points for our scientist customers 
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Example: Dei 


1 


Mission Objective 

Return a large boulder from Deimos 

Launch Vehicle 

Delta IV Heavy with lunar flyby (C3 2.0) 

Power System 


Array power at 1 AU 

chosen by optimizer 

Cell performance model 

1/r 2 

Spacecraft bus power 

2.0 kW 

Power margin 

0% 

Propulsion System 

Thruster 

chosen by optimizer (high-lsp or high-thrust versions of a large Hall 


thruster) 

Number of thrusters 

chosen by optimizer (3, 4, 5, or 6) 

Duty cycle 

90% 

Propellant tank 

unconstrained 

Mission Sequence 

Direct travel to Mars followed by direct return to C3 2.0 for lunar flyby 


capture 


Mars arrival/departure is modeled using an Edelbaum spiral 

Inner-Loop Objective Function 

Maximize sample return mass 

Outer-Loop Objective Functions 

Sample return mass 


Solar array size 


Launch epoch 


Flight time 
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Flight time (years) 


Deimos Samp 
First Generate 




8000 


7000 


6000 


5000 


4000 


3000 


2000 


1000 


0 
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Final journey mass increment (for maximizing sample return) 


Flight time (years) 


Deimos Sampl 
Final Generatk 




10500 


9000 


7500 


6000 


4500 


3000 


1500 


0 


NAVIGATION & MISSION DESIGN BRANCH, CODE 595 

NASA GSFC 


26 * 


Final journey mass increment (for maximizing sample return) 


Deimos Sam; 
Two Trajecto 



LT rendezvous 
Mars 

10/2 S/2021 
m - 8469 kg 


launch 

Earth 

1/1/2019 

C3 = 2.000 km 2 /s 2 
DLA - 18.8 
m = 11191 kg 


departure 

Mars 

6/24/2023 
m = 11900 kg 


le8 


x (km) 

5 - 


begin spiral 
Mars 

10/23/2022 
Ar = 1,079 km/ s 
rn = 12376 kg 


intercept 

Earth 

9/23/2025 

= 1 >414 krn f it 
DEC = 30.5 

m = 10226 kg 


end spiral 
Mars 

5/26/2022 

Av — 1.079 km/s 

m = 8123 kg 


le8 


x (km) 



LT rendezvous 
Mars 

8/18/2021 
m = 9651 kg 


launch 

Earth 

6/18/2020 

C3 = 2.000 km 2 / s 2 

DLA = -9.6 

m = 11191 kg 

end spiral 
Mars 
1/2/2022 
Av = 1.079 km/s 
m = 9304 kg 


departure 

Mars 

12/24/2022 
m = 18974 kg 


intercept 
Earth 
8/31/2025 
r x = 1.414 km/s 
DEC = 25.2° 
m = 16368 kg 

begin spiral 
Mars 

6/18/2022 
At* = 1.079 km/s 
m = 19684 kg 


A 7-year mission with a 41 kW solar 
array returns a 4 ton boulder 


A 5-year mission with a 68 kW solar 
array returns a 10 ton boulder 


Both of these candidate missions, and many others, would be 
of interest to our customers 
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Conclusions 


■ The low-thrust interplanetary mission and systems design problem may be posed as a multi- 
objective hybrid optimal control problem 

■ The combination of a multi-objective discrete NSGA-II outer-loop with a MBH+NLP inner-loop is 
a very powerful way to explore a mission and systems trade space in an efficient, automated 
manner 

■ The algorithm described here has revolutionized the low-thrust interplanetary mission design 
process at NASA Goddard Space Flight Center 

- We can now study multiple mission design cases simultaneously, limited only by available 
computing power 

- Mission design engineers can now spend more time with the customer and with spacecraft 
hardware engineers so that we can fully understand the scientific and engineering context of 
our work 

- Good mission ideas are much less likely to be rejected due to lack of time to work on mission 
design, and bad ideas are much more likely to be rejected before they consume too many 
resources 

■ Skilled analysts are expensive. With a multi-objective HOCP automaton, analysts can focus on 
understanding the customer’s needs and the spacecraft’s capabilities and also detailed design 
work, leaving repetitive tasks to the computer 
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EMTG is available open-source at 
https://sourceforqe.net/proiects/emtq/ 
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ffljC totonic Basin Hopping 
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Tuning Mono 




g (MBH) 


■ We examined two components of classical MBH: 

the random hops are driven by a uniform probability distribution; hops can 
occur in a ball of some user-defined radius about the current best point 

There is a concept of “impatience” - a certain number of iterations where the 
solution does not improve, after which the algorithm resets 

■ In this work we consider: 

Alternative probability distributions (Gaussian, Cauchy, Pareto) which have 
the ability to “hop” over the entire solution space 

Given the above, that the concept of “impatience” may not be necessary when 
using alternative probability distributions 

■ Our objective was to find a version of MBH that would be: 

- Efficient (find better solutions in less time) 

Robust (work well on highly constrained problems and not be sensitive to 
tuning parameters) 
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Probability 

Parameters 


Distribution 

RV Generator 

Excursion 

Parameter 

Uniform 

2 p(r — 0.5) 

p: ball size, 
impatience 

Gaussian 

„2 

5 r 

;=e 2o- 2 

(TV 271 

o\ standard 
deviation 

Cauchy 

p tan(7r(r — 0.5)) 

p: scale 

Pareto 

5 (<x — 1.0) 

6 fchT 

a\ “parameter” 


r = uniform (0.0, 1.0), s is a fair coin flip, e = 1.0 * 10 13 
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Time series of uniform RVs.step = 0.1 



Time series of Cauchy RVs.scale = 0.Q0Z 



Time series of Gaussian RVs, sigma = 0.1 

1.0 


oo iJi Jitfi LliJ jy min litfli i ( iJJ iiL 

■ pWlT^T- f TT| 


-up 

0 200 400 000 800 1000 1200 

Time series of bi-polar Pareto RVs,alpha = 1.01 


1.0 



“ 1.0 


0 200 400 600 800 1000 1200 


We want a distribution which takes lots of small steps to “exploit” the 
local region but also takes frequent large steps to “explore” the rest of 
the space. 
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The Experim 


■ 16 four-day (10000 step) runs of EMTG were conducted for each 
distribution 

■ Each of the 16 runs had a different value of the excursion parameter 

■ Impatience was turned off, i.e. MBH was never allowed to reset during 
the experiment 

- This was necessary to see how effectively each distribution could 
random-walk around the decision space 
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Results - B 
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Results - Av 



-0.35 J 

Uniform (average path, averaged across its 16 excursion parameters) 
Pareto (average path, averaged across its 16 excursion parameters) 
Cauchy (average path, averaged across its 16 excursion parameters) 
Gaussian (average path, averaged across its 16 excursion parameters) 


Pareto-driven MBH is most efficent (better solution in less time) 
and most robust (insensitive to tuning parameters) 
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Why? 


■ Random walks (RWs) can be compared in terms of mean squared 
displacement (MSD) 

- A higher MSD means that an RW travels the problem space faster and 
more thoroughly than a lower MSD 

■ MSD can be used to describe RWs as diffusions through media 

■ In diffusion through homogeneous media (i.e. unconstrained problem 
spaces), RWs driven by independent identically distributed (i.i.d.) 
distributions with finite variance are considered “normally diffusive” 

- MSD proportional to the number of steps 

■ RWs driven by i.i.d. distributions with infinite variance are “super- 
diffusive” 

- MSD proportional to the number of steps raised to some power 
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Why? Continu< 


■ 


■ In a simplified test problem the Pareto RW is super-diffusive while the 
uniform and Gaussian RWs are normally diffusive (or just barely super- 
diffusive for Gaussian) 

■ It is difficult to plot MSD of the Cauchy RW on the same graph because 
Cauchy distributions do not have a mean 

Mean Squared Displacement [step] 
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■ Constraints introduce serial negative auto-correlations 

- Constraint effectively restricts RW from moving in a certain direction, i.e. into the 
constraint 

■ Stochastic global search in constrained problem spaces can be described as diffusions 
through in-homogenous media 

■ When constraints are added to the simplified test problem, the uniform and Gaussian 
RWs become sub-diffusive but the Pareto distribution is still super-diffusive 


Mean Squared Displacement [step] 

based on averaging 200 jndepedent paths through In-htHTioagfious space (with internal <cnstrslfilsi 



— UfiilGriTi, ite-pme - b-4 

— B'-polar Pareto, *lph* n 1 05 
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